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Planetary nebulae as mass tracers of their parent galaxies: 
biases in the estimate of the kinematical quantities 
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Abstract. Multi-object and multi-fiber spectrographs on 4 and 8 meter telescopes make it possible to use extra- 
galactic planetary nebulae (PNe) in the outer halos of early-type galaxies as kinematical tracers, where classical 
techniques based on integrated stellar light fail. Until now, published PNe radial velocity samples are small, with 
only a few tens of radial velocity measurements (except for a few cases like NGC 5128 or M31), causing uncertain- 
ties in the mass and angular momentum estimates based on these data. To quantify these uncertainties, we have 
made equilibrium models for spherical galaxies, with and without dark matter, and via Montecarlo simulations 
we produce radial velocity samples with different sizes. We then apply, to these discrete radial velocity fields, the 
same standard kinematical analysis as it is commonly done with small samples of observed PNe radial velocities. 
By comparison of the inferred quantities with those computed from the analytical model, we test for systematic 
biases and establish a robust procedure to infer the angular momentum distribution and radial velocity dispersion 
profiles from such samples. 
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1. Introduction 

The dynamics of the outer regions of early-type galaxies 
have been studied by means of test particles like globular 
clusters (GCs) and planetary nebulae (PNe). GCs have 
been successfully used to probe the gravitational potential 
in giant ellipticals ( M87: Mould et al. |199C|; C ohen |200C ; 
NGC5128: Sharpies |l98| Harris e t al. |1988| ; NGC1399: 
Gri llmai r et al. |1994| ; Minniti et al. |1998 ; Kissler-Patig et 
1998| ), but their number density and angular momen- 



al 

turn distributions often turn out to be different from those 
of the stellar population in the outer galaxy halos. 

PNe are a population of dying stars, whose outer en- 
velope re-emits more than 15% of the energy emitted in 
the UV by the internal star in the [OIII] green line at 5007 
A (Dopita et al. 1992) and therefore they can be readily 
detected in distant galaxies. The observational evidences 
indicate that the number density of PNe, unlike GCs, is 
proportional to the und erlyin g stellar light (Ciardullo et 
al. |1989|; C iardullo et al. |l99l| ; McMillan et al. |l993| ; Ford 
et al. 1996|) and they share the angu lar m omentum distri- 
bution of the st ars (Arnabo ldi et al. 1994: Hui et al. 1995 ; 



Arnaboldi et al. 1996, 1998). First attempts to use PNe as 



to 1986 (Nolthcnius & Ford |1987| ) and studies of early- 
type systems within a distance of 10 Mpc rapidly followed 
(Ciardullo et al. Tremblay et al. Hui et al. 



1995). Sinc e 1993 new observing techniques (eg Arnaboldi 



et al. 1994) allowed measurements of radial velocities of 
PNe in the outer regions of giant early-type galaxies sit- 
uated at distances larger than 10 Mpc. These studies, 
mostly based on samples of only a few tens of PNe, show 
that the outer PNe typically have faster systemic rotation 
than GCs (Gri llmair et al. [l994 Arn aboldi et al. |1994 ; 



Hui et al. 1995: Arnaboldi et al. 1 



kinematical tracers in nearby galactic systems date back 



So far these small samples were analysed by adopt- 
ing simple three-parameter functions for the underlying 
projected rotation field, but no detailed studies have yet 
been reported to test for biases in the estimated rota- 
tion and velocity dispersion introduced by these adopted 
(parametric) functions. Non-parametric analyses can in 
principle avoid this problem. In practice, however, for the 
small data samples of interest here, the associated inher- 
ent smoothing effectively drives the non-parametric anal- 
ysis towards one of the commonly used simple parametric 
forms (like solid body rotation: see for example Arnaboldi 
et al |199^). 
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In the large telescope era, much larger samples of PNe 
velocities will become available for galaxies out to dis- 
tances of about 20 Mpc. The question of biases induced 
by adopting simple parametric forms for the rotation field 
remains relevant, because studies of galaxies at larger dis- 
tances will still be limited to small sample of radial veloc- 
ities. It will be important to know how to compare the re- 
sults for nearby giant ellipticals, based on large samples of 
PNe velocities, with those for more distant objects, based 
on smaller samples. 

In this paper we build simple equilibrium models with 
and without dark matter (DM) for spherical early-type 
galaxies (Sect. Via Montecarlo simulations, we extract 
samples of tracers of various sizes (Sect. |^). We then adopt 
simple forms for the mean rotation field and derive the pa- 
rameters for the projected 2-D rotation velocity field from 
these samples (Sect. 4 and 5). We also derive estimates of 
the precision and biases for the mean rotation and velocity 
dispersion (Sect. ^ and 7). The main goal of this work is to 
determine the minimum sample needed to derive reliable 
estimates of the kinematics of the host system via these 
parametric fits, and to define a robust approach to derive 
observables from small radial velocity samples. Discussion 
and conclusions arc drawn in Sect. 
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Fig. 1. Mass distributions for the Hernquist model: luminous 
mass matter (dash-dotted line), dark halo mass (dashed line) 
and the total mass (solid line). 



1.114. Writing r for r/a, we define a dimensionless density 
distribution 
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2. Analytical spherical systems in equilibrium 

In our modeling we wish to represent the characteristics 
of the observed PNe distribution in early- type galaxies. 
In these systems, the observed PNe number densit y fol 
lows the stellar light distribution (CiarduUo et al. 19i 
CiarduUo et al. |l99l| ; McMillan et al. |1993| ; Ford et al. 
1996), except in the bright inner regions where the PNe 
counts become incomplete. As the surface brightness of 
the stellar continuum background increases towards small 
radii, the 5007 A [OIII] emission from the PNe is more 
difficult to detect and the apparent ratio of PNe to lu- 
minosity decreases. In Appendix, we compute the value 
of the limiting radius i?iim for which at large radii the 
PNe sample is complete. For the surface brightness pro- 
file of a typical E galaxy in Virgo Cluster, we find that 
-Riim = 0.7 X i?e, where i?e is the effective radius. 

2.1. Systems without dark matter 

Assuming const ant m ass to light ratio (M/L), the ana- 
lytic Hernquist ( 1990| ) model is a good approximation to 
a system whose surface brightness distribution follows the 
de Vaucouleurs law ( 194§| ). The luminous mass density is 
given by: 



p(r) = C\ 



M\a 



27r r[r -t- aY ' 



(1) 



where M\ is the total luminous mass, a is a distance scale 
(i?e = 1.8153 a) and C\ is a normalization constant. We 
consider systems truncated at i?max = 18a (i.e. ~ 10i?c). 
For a = 1 and M\ = 1, the normalization constant C\ — 



The cumulative mass distribution is then 

M\{r) = An I x^p{x)dx. 
Jo 



(3) 



The test particles are extracted from this distribution, 
taking into account the incompleteness effects at radii 
< i?iim- We then consider the kinematics of PNe within 
10a = 5.5i?e, corresponding to the typical radial coverage 
of the observed PNe samples in real galaxies. 

2.1.1. Systems with dark matter 

In constructing our equilibrium models, we also consider 
systems with an additional mass contribution to the stellar 
mass density coming from a dark halo. Furthermore we 
model the dark halo with a Hernquist mass distribution, 
with scale length large enough so that it roughly mimics 
a p ~ r^^ halo in the region of interest (r < 10). We write 
the dark matter density as 



Pdir) = 



1 



27r r(r -|- d)^ 



(4) 



and adopt d = 10a and — 7 .1M\ (in agreement with 
the mass distribution of NGC5128, Hui et al. |199 



cumulative mass distribution is defined as in Eq. 
the potential is derived from the total mass 

Mt(r) = Afi(r)+Md(r). 



The 
and 



(5) 



The Mi(r), Md{r) and A/t(r) distribution are shown in 
Fig. I 
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2.2. Intrinsic and projected kinematics 

We consider non-rotating and rotating equilibrium spheri- 
cal systems, with the velocity dispersion components given 
by the Jeans equations^]. 

For non-rotating systems with no dark matter, we assume 
that mass follows light and adopt an isotropic velocity 
dispersion to solve the radial Jeans equation 



dr 



-G 



Mi{r)pir) 



(6) 



The solution in this simple case is 
cr2(r) = 7.110-^r(l. + r)3{-0.13- 



[130321(1 



^25 + 52r + 42r2 + 12r3 (7) 



+ 12(1 + r)4 log(r) - 12(1 + r)^ log(l + r))]} 
where a{r) is given in units of Vg = VGMi/af|). 

We should also investigate rotating systems because 
early-type galaxies appear to show fast rotation in their 
outer parts (Hui et al. |1995| ; Arnaboldi et al. |1996| , |1998D . 
What is the most appropriate form for the rotation law 
? The only system for which a mean stellar rotation law 
has been derived re liably out to large radii is Centaurus 
A. Here, Hui et al. ( 1995 ) adopt the functional form 



l'rot(^) 



(8) 



where R is the distance from the galactic center and i?v a 
scale distance, which reproduces the rotation of NGC5128 
in its equatorial plane. We assume that our model systems 
have the same underlying rotation law, and adopt a cylin- 
drical rotational structure with the same functional form, 
where R is now the radial coordinate in a cylindrical (i?, 
V?, z) system. 

In a conventional cylindrical coordinate system, (i?, (p, 
z), the Jeans equations become 



OR 



dpvl 9$ 

=0, 

az oz 




and 



(9) 



(10) 



where p and $ are functions of (i?, z). 

If we substitute into equation (H) a given cylindrical 
rotation law, and assume an isotropic velocity dispersion, 
then this system is over-determined. If we now relax the 
condition of isotropy and assume that the principal axes 



^ The possibility that s pheri cal systems can rotate was inves- 
tigated by Lynden-Bell (196C): here we consider such systems 



in order to use simple density and gravitational potential func- 
tions. 

^ If we take Mi = 3.5 x lO^^Af© and a = 2 kpc, we have 
Vg — 870 km/s, so the velocities in our simulations are of the 
order of 0.3«g — 260 km/s. 



of the velocity ellipsoid lie parallel to the (i?, z) axes and 
(7^ = crj ^ cr^ ; the System of equations can be written as 



dp{R,z)al(B.,z) 
dR 



p{R., 



dR 



dp{R,z)al{R,z) 
dz 



= -p{R,z) 



d^{R,z) 
dz 



(11) 



(12) 



These two equations can be solved independently by sim- 
ple inversion as in Binney & Tremaine (1987)^. We believe 
this to be a reasonable approach for our purpose. It de- 
scribes an equilibrium system with a spherical potential 
and a given cylindrical rotation law of known functional 
formQ. 

Physically plausible solutions require that the compo- 
nents of the velocity ellipsoid are all positive in the volume 
occupied by the system. As we will see later, only certain 
values of V^nax in Eq- ^ allow physical solutions in the 
whole volume occupied by the system. In Fig. || we show 
the intrinsic kinematics for Vmax = 0.23wg and i?v — 2a. 

2.3. Geometry and projected kinematics 

As the galaxy intrinsic coordinate system we adopt a 
Cartesian coordinate system where the X-axis coincides 
with the direction of the line-of-sight and the Y-Z plane 
is the sky plane; the positive direction of the X-axis is to- 
wards the observer, while VJ-ad — —Vx is the radial veloc- 
ity. In case of inclined systems, we adopt the conventional 
Euler matrix to identify the coordinate system X'Y'Z', 
where the X'-axis is coincident with line-of-sight and the 
Y'-Z' plane is the Sky plane. The three Euler angles are 
0, e, and (Goldstein, H. \198(\): when ,^ = 6* = -0 = 0, 



i.e. the edge-on case, the two coordinate systems XYZ and 
X'Y'Z' are coincident. 

The projected rotation and velocity dispersion profiles are 
computed from the intrinsic profiles as follow. The pro- 
jected moments of the radial velocities along the line-of- 
sight are given by the Abel integral 



(^os ) Tf^\ I ^ad / 



pr 



-dr 



(13) 



where t is the distance from the galactic center along 
a fixed direction on the sky plane, I(t) is the projected 
density, and is the n-th moment of the radial veloc- 
ity in each volume element along the line-of-sight. When 



^ They solve the second equation in the isotropic case (page 
120). 

^ Our assumptions on the alignment and anisotropy of the 
velocity dispersion tensor were motivated by the need for a 
simple approach to the solution of the Jeans equat ions. Our 
choice is also physically motivated by Arnold (1995) who ob- 



tained an axisymmetric solution which is quasi-cylindrically 
aligned and where a\ ~ ff^ 7^ (model vi). An advantage 
of this approach is that an analytical solution can be obtained 
and used during our simulations. 
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Fig. 2. R and z dependence of the velocity dispersion ellipsoid (left) and anisotropy parameter P — 1 — cTz/'^r (right) for the 
rotating model with DM and V^ax = 0.23wg. In the left panel the contours of the ctb(= a^) component ((Tj?/iig=0.325 to 0.125 
with step=0.025 from inside) are plotted with different grey scales; the circular lines are the corresponding contours. 



we consider inclined systems, we have V^ad ~ ^^'^ = 
{aiiVx + ai2Vy + ai^Vz) in Eq. { [T^ ) where are the ele- 
ments of the Euler matrix. Vrad can be expressed in terms 
of the cylindrical velocity components: 



Vx = Wi?cos</? — VtpSmLp 
Vy = vr sin Lp -\- cos Lp 



(14) 



where ip is the conventional azimuthal angle in spher- 
ical and cylindrical coordinate systems. Then, VJ-ad = 
-{ciVR + C2V^ + czVz) where 



Ci = flu cos + sm 
C2 — — flu sin -f ai2 cos 
C2 = ai3- 



(15) 



When Eqs. ( |14D and ( |15D are substituted in Eq. (|13|), we 
derive the projected kinematics along the line-of-sight. 
The rotation curve expected from the model at equilib- 
rium is given by (ViosCi")), and the velocity dispersion pro- 
file is given by 4 ^ (y,l,{T)) - (Mos(r))2. 

3. The simulated discrete radial velocity fields 
(SDRVF) 

The dynamical state of a galactic system is described by its 
phase space distribution function (DP). With a functional 
form for the density distribution and gravitational poten- 
tial, it is possible to determine the DP via the Eddington 
formula (Binney & Tremaine 1987) for spherical systems. 

In the case of the adopted Hernquist model, the 
Eddington formula is not so useful because its applica- 
tion required the use of special functions. To deal with 
the DP in a simple heuristic way, we adopt a factoriscd 
form 



where p(r) is the mass density distribution given by the 
model and -F(r, v) is the velocity distribution depending 
on the position, r, via its moments. 

The observed distribution of radial velocity profiles in 
early-type galaxies are Gaussian, with maximum deviation 

Bender et al. 



of 10% (Winsall & Preeman 1993 



1994). In 



our modeling we will assume a velocity distribution given 
by the product of three Gaussian functions, although some 
anisotropy is implied by the assumed intrinsic cylindrical 
rotation. We write the velocity distribution as 



where 



Fi{Y,Vi) 



2Tl(Ji{v) 



exp 



(«.-V-(r))2 



(17) 



i = z(18) 



and Fr, F^p, F^ are the normal distributions along the 
directions of the principal axes of the velocity ellipsoid. 
Por each Gaussian, the mean value Vi represents the mean 
motion and the standard deviation ai is the diagonal 
element of the velocity ellipsoid in the same direction, 
both evaluated at position r 

We consider both non- rotating and rotating systems: 
once the rotation velocity structure is assigned, we solve 
the Jeans equations for the given density distribution and 
obtain the velocity dispersion ellipsoid, plus the related 
projected velocity dispersion profiles. In rotating systems, 
the velocity distribution is given by Eq. (^, setting 
Vi = v-cot as in Eq. fH) and ai as solutions of the Jeans 
equations Eqs.(|ll|), (|l^); in non-rotating system, Vi — 



/(r,v) = /9(r)F(r,v) 



(16) 



^ The DF as given in Eq. (^, with the assumptions (Q 
and ( p^ for the velocity distribution, is not a solution of the 
coUisionless Boltzmann equation and does not represent the 
true dynamical state of a stellar system. However, for our pur- 
pose this approximate representation of the DF is taken as a 
reasonable description of an instantaneous state for a spherical 
system. 
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and ai is given by Eq. 

We now generate a sample of discrete test particle from 
the mass density distribution p{r) and the distribution 
function i^(r, v) as follows: 

1. we randomly extract from p(r) the spherical coordi- 
nates (r, 6, If) for each star and check that the corre- 
sponding projected radius R is greater than i?iim (see 
Sect. 2), 

2. from the velocity distribution F{y, v) we randomly ex- 
tract the three velocity components (wfl, w^, fz) at the 
position r of each star. 

Once the 3D velocity field is obtained, we project it 
on the sky plane. The next step is to simulate a measure- 
ment of this two-dimensional velocity field: 1) we assume 
a typical measurement error for low to medium dispersion 
multi-object spectroscopy (cTmca ^ 70 km/s — O.OSwg, see 
footnote 1), with a normal distribution; 2) the observed 
radial velocity measurement is obtained by a random ex- 
traction from a Gaussian, whose average is the Wj-ad value 
obtained for a single PN and standard deviation CTmoa- 
In this way we obtain a simulated discrete radial velocity 
field, Vohs{X, Y), from which the kinematical information 
is extracted. 

4. Standard analysis procedure 

From our simulated systems we want to estimate the kine- 
matical quantities without using our knowledge of the in- 
trinsic dynamics of our system. Therefore we apply to the 
simulated data a procedure similar to those adopted in 
the analysis of the small samples of discrete radi al ve- 
locity fields of observed PNe (eg Arnaboldi et al. 1994 



1996| , |1998[ ). In this procedure, we fit some simple three- 
parameter functions to our simulated data Vobs {X^ Y) . We 
do this to see how misleading these fits can be if the real 



We perform a 



rotation fields have the form of Eq. 
last-squares fit to the Yo\ys{X,Y) using: 

1. a bilinear Junction (hereafter BE) 

Wrad = a + bX + cY (19) 

where X, Y are Cartesian coordinates of the PNe (typ- 
ically in CCD coordinates) and a, 6, c are parameters 
to be determined by fitting Eq. ( p^ to the data. These 
parameters are used to determine the systemic veloc- 
ity, direction of the axis of maximum velocity gradient, 
Z\, and modulus of the velocity gradient. This linear 
velocity field is equivalent to a solid body rotation of 
the form Wrad = Wsys + cos{4> ~ 4'*) where Vsys is the 
systemic velocity, uj is the angular velocity, (j) is the 
angle from the kinematic major axis, and 0*, the P. A. 
of the Z\ axis. 

2. a fiat rotation curve (hereafter EC) 

Vr&d = Vsys + COS(0 - (/)*) (20) 

where Vsys and (f)^ are the same quantities as in the 
BE, and is the amplitude of the flat rotation curve. 



The residual field, Aw, is computed as the difference be- 
tween Vobs and the interpolated radial fields from the fits. 
The residual fields obviously depend on the adopted fit 
(BE or EC). Eurthermore we assume point-symmetry to 
the galaxy center for the simulated samples. Therefore 
to every velocity v at position (X, Y) on the sky plane, 
there is a corresponding velocity —v at position {—X, — Y). 
Hereafter we will refer to a symmetrized velocity field when 
we consider a sample generated on one side of the galaxy 
by adding the symmetric points generated by test particles 
on the other side. 



5. Simulations 

In Table |l| we summarise the set of parameters used in 
our models to generate the SDRVEs. Models with V^iax 
up to 0.23z;g are physical, i.e. they have velocity disper- 
sion ellipsoid components which are positive everywhere, 
while for T^nax > 0.23wg the (tu = components become 
negative for large z. However, we have kept these models 
because they have physical solutions near the equatorial 
plane, where the investigation of the kinematical observ- 
ables (rotation velocity and velocity dispersion profiles) is 
often focussed. 

The solutions of the Jeans equations are shown in 
Fig. ^ for the system with l^iax — 0.23wg and for the 
non-rotating system. For our models we produce samples 
with 50, 150 and 500 PNe. For a given sample size, a fixed 
mass density model, a particular set of values for the ro- 
tation velocity (as given in Table 1) and inclination of the 
system, we produce 100 realisations of the radial velocity 
fields as described in Sect. 3. Fig. ^ shows a view of the 
the radial velocity fleld and the projected number density 
distribution obtained for a sample of 500 PNe for a ro- 
tating simulated system with Vmax = 0.23tig (see Table |l]) 
seen edge-on. The comparison with the luminous matter 
projected density distribution shows the expected incom- 
pleteness effect in the generated PNe sample towards the 
center of the mass distribution. When computing the ex- 



pected model quantities from Eq. (13), we take this selec- 



tion effect into account, using the spatial and projected 
density profile obtained from a sample of 50 000 objects. 



5.1. Global kinematical quantities 

To each SDRVE, we apply the analysis procedure de- 
scribed in Sect. ^ and derive the following kinematical 
quantities: 



1. the velocity gradient, gradV = ±V&^ + for the BE, 
with sign given by the sense of rotation; 

2. the maximum rotation velocity (the parameter) for 
the EC; 

3. the position angle of the maximum gradient, ^zi ~ 
arctan(c/6) for the BF, or ^zi = 0* for the EC. 
Hereafter we refer to the axis of the maximum velocity 
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Fig. 3. Here we show the solutions of the Jeans equations (Eqs. and |l2|) as functions of distance from the center (R=r, z). 
Left: isotropic non-rotating system without DM. Right: rotating equilibrium model with DM and Knax = 0.23115 (see Table 0); 
the o"R = ay, component is drawn as solid line on the equatorial plane and as long-dashed line along the rotation axis, while the 
(Tz component is the dash-dotted line. The rotation velocity on the equatorial plane is also plotted, as short-dashed line. 




Fig. 4. Left: a typical SDRVF obtained for a sample size of 500 PNe and the equilibrium model with DM and cylindrical 
rotation (inclination angles: 0:0:0). The dashed line is the line of the maximum gradient (computed with the bilinear fit). 
Right: the PNe number density (histogram) is compared with the surface brightness of the simulated system (solid line). The 
incompleteness of the sample is evident in the inner regions. 



gradient as Zl, with Z2 its perpendicular axis, i.e. the 
apparent rotation axis; 



4. the systemic velocity Vgys 
rameter from FC. 



a for the BF and 



pa- 



For each set of 100 simulations and a given sample size, 
we obtain the distribution of these global quantities, the 
related mean value and standard deviation of the sample 
(SD)^, which are then compared with the expected values 
to check for any presence of biases. 



In this work we will use the abbreviation SD instead of the 
usual symbol a in order to avoid confusion with the velocity 
dispersion. 



5.2. Projected kinematical quantities 

To investigate the properties of our SDRVFs, we will im- 
plement a binning procedure, because similar strategies 
were applied to small observed radial velocity samples (like 
those for NGC 1399 and NGC 1316) and also in the case 
of NGC 5128 for which 433 radial velocities were available. 
Similar procedures are still adopted when using GCs data 
(eg Minniti et al. [1998|; Kissler-Patig et al. |1998D. In this 



framework we can also establish whether there are biases 
introduced by binning the observed radial velocity fields. 

Once the line of maximum gradient is identified, the 
symmetrized velocity field is spatially binned along Zl, se- 
lecting particles within a strip along Zl. We have adopted 
strips having different width dZ in order to check whether 
there was some dependence of the kinematical estimates 
on this parameter: see Table 3. The spatial binning is done 
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Table 1. Parameter values for the dynamical models 
which are used to generate the 2-D discrete radial veloc- 
ity samples. The a and Mi parameters are common to the 
self-consistent non-rotating and dark matter rotating sys- 
tems. Here we adopt different values for the rotation vs. 
velocity dispersion so that any bias effects caused by the 
assumed functional form of the rotation field can be in- 
vestigated. We refer to Sect. ^ for a detailed discussion. In 
parentheses we show the v/a value at i? = 10a. 



Dynamical and kinematical parameters 


a 


1 


Ml 


1 


d 


lOa 




7.7Mi 




2a 




0.35 {v/a^l.8) 




0.30 Dg (1.4) 


Vmax 


0.23 Dg (1.0) 




0.20 Wg (0.8) 




0.15 Vg (0.6) 



in such a way that the number of particles in each bin is 
about 10 or more. We also do a radial binning by selecting 
PNe in radial annuli. This binning is usually done for ra- 
dial velocity samples which are small, assuming spherical 
symmetry and isotropy[]. In our analysis of a given sample 
size, the radial and Zl binning are fixed and maintained 
for all the 100 SDRVFs (the same binning is adopted for 
the Z2 axis because of spherical symmetry). 

In the standard procedure, the mean values of the 
radial velocity sample in spatial bins along Zl provides 
a measure of the rotation curve Vrot, while the velocity 
dispersion is obtained as SD in the same bins from the 
residual fields for the BF or FC. This gives an estimate 
of the observed dispersion (Tobs from which the projected 
velocity dispersion is obtained, using the independently 
determined measuring error amca, as try = V^'^bs^^'mca- 
This velocity dispersion obviously depends on the adopted 
functional form of the rotation field. For comparison, we 
also estimated the velocity dispersion independently of the 
adopted rotation field, as the SD of the velocities from 
the mean radial velocity in each bin. Hereafter we refer to 
these estimates as NFP (no- fit procedure). We associate 
to each value of mean radial velocity Vrot an error given by 
^^rot = SD/\/N and to dobs an error Suohs — SD/\/2N . 
The behavior of Vohs, cry and their relative errors from 
100 simulations allow us to estimate the precision of any 
observable for a given sample size. 



6. Results 



Here, the analysis on radial bins is performed only for mod- 
els which are everywhere physical (Knax < 0.23wg). 



6.1. Estimates of Kinematical Parameters 

In Table || we list the results from our fits for the rotating 
model seen from three different line-of-sights (as identified 
by the Euler angles). For each set of simulations and fits. 
Table ^ shows the expected and estimated values of the 
kinematical parameters. The systemic velocity, Vsys, and 
the P. A. of the maximum velocity gradient, $21 , are well 
estimated from the fits of a simple three-parameter func- 
tion. They are consistent with the expected model values 
at the 95% confidence level, for the potential with and 
without dark matter, and for both functions (BF and FC) 
that were fit to the velocity fields. We note the impor- 
tance of an accurate estimate for the position angle 
because ^zi is then used in the analysis of the projected 
kinematics. Systematic errors on ^zi produce systematic 
errors on the projected kinematics. 

The velocity gradient and maximum rotation velocity 
estimates are systematically incorrect when they are ob- 
tained from the wrong simple function. In the case of cylin- 
drical rotation, the BF can produce a large over-estimate 
of the velocity gradient (up to about 100%) while the FC 
(which is closer to being the correct functional form for the 
rotation adopted in the model) gives only a small under- 
estimate. The reason is clear from examination of Fig. ^ 

6.2. Estimates of l<inematical quantities 

Rotation velocity curves The mean velocity rotation 
curves evaluated from 100 SDRVFs for each given sample 
size follow the theoretical behaviour predicted by the 
models for every parameter choice (Fig. |l^). Rotation 
along Z2 is everywhere consistent with the null rotation 
(Fig. |l^ 

In Figs nq and y we show the results obtained considering 
the inclined system: the rotation velocity is decreased by 
the expected sini factor, where i is the angle between the 
line-of-sight and rotation axis given by a particular set of 
Euler angles. For face-on model, no-rotation is found, as 
expected (Fig. p] ). 

Velocity dispersion profiles along the Zl and Z2 axes - 
The inferred velocity dispersion profiles are of particular 
interest because they allow the evaluation of the mass 
distribution via inversion of the Jeans equations. When 
the velocity dispersion profiles are computed from the 
residual velocity fields, they can be affected by biases 
resulting from the simple functional forms used to fit the 
discrete radial velocity sample. 

When the equilibrium model with cylindrical rotation 
is compared with the results for larger samples (eg 150 
and 500), the velocity dispersion profiles obtained from 
the BF residuals along the Zl axis show a small but vis- 
ible over-estimate of cry in the outer bins relative to the 
expected values from the model. Better agreement is found 
using the FC (see Fig. This happens because, along 
the maximum gradient axis, the flat rotation curve gives a 
better description for the rotational structure of the equi- 
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Zl/a Zl/i 



Fig. 5. Biases in the maximum gradient and maximum rotation velocity. The intrinsic rotation curve from the model (curved 
solid line, Knax = 0.3ug, see Table ^ is compared with the bilinear fit (left), the flat-curve fit (right) and the mean values in 
bins (flUed points) obtained on 100 simulations. The parameters of the fitted functions are the averages from the same 100 
simulations: straight solid lines show the BF (left) and FC (right) fits, and dot-dashed lines are the ±SD lines about the fits. 
Quoted errors about the points are the SDs of the distribution of estimates in each bin. The dashed vertical line indicates the 
inner limiting radius, i?iim- 

Table 2. Results obtained for the equilibrium model with DM and flat cylindrical rotation: behavior from 100 SDRVFs. 
Square parentheses: the unit in which the quantities are given; round parentheses: the expected values. Quoted errors 
are SDs of the derived quantities. 



Model with DM - Rotation Curve of Equation (9) 



Proced. 




BF 






FC 




Inclin. 


</> = e = ?/> = 


N° PNe 


<&zi[rad](0) 


gradV[ug/a] (0.023) 


Vsys[v,] (0) 


$zi[rad](0) 


v4v^] (0.23) 


V.y.K] (0) 


500 


0.00±0.08 


0.046±0.005 


0.002±0.013 


-0.01±0.09 


0.182±0.019 


0.000±0.014 


150 


-0.06±0.22 


0.045±0.010 


0.01±0.02 


0.00±0.23 


0.18±0.04 


0.00±0.02 


50 


0.0±0.3 


0.048±0.015 


0.0±0.04 


0.0±0.3 


0.18±0.06 


0.00±0.04 


Inclin. 


(j) = -45 6 = -45 V = -45 


N° PNe 


$zi[rad] (0.62) 


gradV[i;g/a] (0.021) 


V^ysK] (0) 


$zi[rad](0.62) 


v4vja] (0.21) 


VsysK] (0) 


500 


0.61±0.12 


0.039±0.005 


-0.002±0.013 


0.62±0.12 


0.155±0.022 


-0.003±0.014 


150 


0.61±0.23 


0.040±0.009 


0.00±0.02 


0.61±0.23 


0.16±0.04 


-0.01±0.02 


50 


0.5±0.4 


0.041±0.013 


0.00±0.04 


0.5±0.3 


0.16±0.05 


0.00±0.04 


Inclin. 


<;!) = 90 6* = 90 i/- = 30 


iV° PNe 


<l>zi[rad] (indet.) 


gradVbg/a] (0) 


VsysK] (0) 


"I>zi[rad] (indet.) 


v4vja] (0) 


VsysK] (0) 


500 


0.13±0.99 


-0.001±0.006 


-0.002±0.014 


0.14±0.67 


O.OliO.Ol 


-0.002±0.014 


150 


0.01±0.89 


0.001±0.013 


0.00±0.03 


0.06±0.62 


0.02±0.03 


0.00±0.03 


50 


0.09± 0.98 


0.000±0.022 


0.00±0.05 


0.18±0.65 


0.04±0.04 


0.00±0.05 



librium model under study, whose rotation becomes flat 
at large distances from the center along the Zl-direction. 

We flnd a clear and expected correlation of this bias 
with the v/a ratio. In the last bin along Zl, we com- 
pute the difference Act = ay — CTmod, between the dis- 
persion values estimated from BF (ctv), and the expected 
dispersion from the models (CTmod), for different values of 
(w/CT)obs = K-ot/cv values in the last bin (obtained by 
varying the K parameter in the rotating models). The 
results for 500 PNe are shown in Fig. ^. The bias in the 
velocity dispersion in the last bins is an increasing function 
of the (f /CT)obs ratio: larger rotation causes the bilinear flt 
to introduce a larger overestimate of the velocity disper- 
sion proflle. 

As mentioned earlier, for small radial velocity sam- 
ples this kind of bias is associated not only with our 
simple parametric flts but also with more sophisticated 



non-parametric analyses. When there are fewer data 
points, the non-parametric algorithms with their inher- 
ent smoothing effectively flt a plane through the data, so 
the velocity dispersion proflles derived from the residual 
flelds do suffer from similar biases as the one we are in- 
vestigating in the case of a simple bilinear flt. 

In general, for all sample sizes, the velocity dispersion 
estimates along the Zl axis obtained using the NFP are in 
better agreement with the expected ones. The estimates 
along Z2 axis are not affected by this bias: in Fig. |l^ 
the estimates from BF and FC are in agreement with 
the model, while NFP systematically overestimates the 
velocity dispersion values in the bins. This effect depends 
on the strip dimension adopted to select particles along 
the rotational axis: the expected proflle is computed along 
the Z2-axis, while the estimates are obtained on bins and 
their widths include regions where the velocity dispersion 
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7. Sampling errors 




(v 

Fig. 6. Biases from BF in dark matter systems: points 
represent the relative deviation of the velocity dispersion 
determined in the last bin along Zl with respect the "true" 
value expected from model. Error bars take in to account 
errors on the estimates on 100 simulations. The behavior 
shown by the residuals indicates that higher rotation in 
the model causes the bilinear fit to introduce a larger bias 
in the estimated velocity dispersion profile. The vertical 
arrow shows the bias in the velocity dispersion estimate, 
for the model with Vniax = O.SOwg (see also Figure 



is larger (see Fig. g). On average, they trace a larger 
dispersion and this causes an overestimate. On the other 
hand, the different fit procedures subtract a gradient in 
these bins, and produce a lower dispersion. In Fig. |lj, the 
estimates with the smallest dZ2, adopted for the different 
samples, show the magnitude of this effect. Smaller dZ2 
reduce the overestimate but the errors are larger. All 
these features are found for the inclined systems too, as 
shown in Fig. Here the profiles follow those of the 
edge-on case, but the amount of bias due to the fitting 
procedure is smaller because the (w/cr)obs = Kot/fv ratio 
is lower^. In the face-on case, no biases are found as in 
case of non-rotating self consistent model (Fig. p^ . 

Velocity dispersion profiles from radial bins - When 
we derive the velocity dispersion profiles from the residual 
velocity field (either from BF or FC) in radial annuli, 
they show effects from the azimuthal dependence of 
the intrinsic velocity dispersion ellipsoid. Those profiles 
derived either from BF or FC have a behavior which 
is intermediate between those expected from the model 
along Zl and Z2 (Figs. |l^ and ^6|) for the edge- on and 
the inclined cases. In the face-on case (Fig. |l^, no 
azimuthal dependence is expected or seen (see Fig. y). On 
the other hand, the velocity dispersion profiles obtained 
with radial binning and NFP show a strong overestimate 
of the velocity dispersion values due to the cylindrical 
structure of the intrinsic velocity field. 



The sampling errors are related to the size of a given sam- 
ple, and affect the accuracy of the kinematical quantities 
we are trying to estimate. They are important because 
they often dominate the error budget of the estimated 
mass and angular momentum estimates in halos of galax- 
ies. In this work, for given sample sizes, we simulated 100 
realisations of the radial velocity field, from which we de- 
rived the distribution of the kinematical measurements. 
We then compare the uncertainty estimated (internally) 
from single measurements of a kinematical observable with 
the uncertainty from the distribution of that observable 
from 100 simulations. This allows us to check i) whether 
the statistical errors are evaluated in a realistic way for 
single measurements and ii) study the behavior of the 
sampling errors and the precisions based on their actual 
distributions. 



Table 3. Average relative precision considering all bins 
and models: RB are radial bins, d.^min and dZj^ax are the 
minimum and maximum width of the strips adopted along 
the Z-axes. In parentheses we indicate the mean number 
of PNe selected with our choice of the bin dimensions. 



Relative 


errors on velocity dispersion 




N° PNe 


RB (iV) 


d^min (N) 




500 
150 
50 


0.106±0.001 (55) 
0.188±0.005 (20) 
0.238±0.016 (15) 


0.211±0.015 (23) 
0.310±0.015 (9) 
0.38±0.03 (8) 


0.164±0.006 (34) 
0.233±0.007 (14) 
0.251±0.015 (12) 



Table 4. Average relative precision in the last bin: RB 
are radial bins, dZ^i^ and dZ^^a.^ are the minimum and 
maximum width of the strips adopted along the Z-axes. In 
parentheses we indicate the mean number of PNe selected 
with our choice of the bin dimensions. 

Relative errors on velocity dispersion in last bin 



N° PNe RB (iV) 



dZmiu (N) 



500 
150 
50 



0.12±0.02 
0.21±0.02 (17) 
0.28±0.02 (11) 



0.32±0.03 (8.5) 
0.37±0.04 (8) 
0.47±0.08 (7) 



0.24±0.03 (12.5) 
0.25±0.03 (13) 
0.29±0.03 (11) 



* The net effect of the inclination is to decrease the 
{v/a)ohs ~ Kot/o"v ratio in the simulated systems. 



7.1. Errors on the velocity dispersion 

For all models and velocity dispersion profiles, we compute 
the average of the relative errors on 100 simulations in all 
bins (for which the PNe number density is complete) for 
a given sample size. We consider this to be indicative of 
the precision for the velocity dispersion. We do this for 
the spatial binning along Zl, Z2, and for radial bins. The 
results are shown in Table pl 
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N (v/tr)... 

Fig. 7. Left: behavior of the relative errors on the velocity dispersion with respect the mean number density in the last bin (full 
points) and in two more bins (see footnote 9) on all the simulation done for all the models. Solid line is a fit to the data points 

using the function A/ViV + B. Right; correlations between relative errors for the rotation velocity and {v/a)obs in the last 
bin. Here difi'erent symbols are adopted for different models, sample sizes and bin dimensions. They fill in the region identified 
by the fitted function for the maximum mean sample, A''=13, and the minimum one, N=7, per bin used in our simulations. 
Dash-dotted lines are for intermediate number samples. 



The average relative errors depend only on sample size 
and are quite independent on the fitting function, inclina- 
tion, and the mass and rotation models^. If we consider 
the estimates along the Z-axes, for a sample size of 500 
test particles the precision on velocity dispersion is 16%, 
for 150 test particles it is 23%, and for a sample of 50 test 
particles it is 25%. Higher precisions (10% for 500 PNe, 
19% for 150 PNe, 24% for 50 PNe) are found if we consider 
the estimates in radial bins (where the whole sample size 
is used). The mean number of particles related to these 
estimates are shown in Table ^ 

The results for the last bins arc shown in Table ^. The 
relative errors are plotted in Fig. Q against lVp"| showing a 
pseudo-Gaussian behavior. The distribution is reproduced 
by the function A/\^ + B, with best fit parameters, 
A=0.94 ± 0.06 and B=-0.01 ± 0.02. If Say /ay = l/Vm, 
the B and A values should be zero and 1/V2 = 0.71 re- 
spectively. Our estimate for B is consistent with zero. We 
find A=0.94>0.71 because by definition Sap/ dp = 1/V2N 

4 " " 



Say 
ay 



^=a^>^ = ^. (21) 
l^p CTp CTp cTp V27V 



For a typical (Tinoa/o'p=0.4 in last bins, Ccr = 1.2 
and the expected value for the A parameter is 
A=C^/^/2=l. 2x0.71=0.85. Furthermore ay = (ct^ - 
^mea)^^^ Is uot a Gaussiau variable, as we do not have 
a linear relation between ay and dp, so this causes an 



The variation of the precisions between models is less than 
10% for a given sample size. 

Here are plotted the results from Tableland additional two 
values obtained for the 150 PNe in the 4"* bin and for 500 PNe 
in the 6*''bin. These values for N=14 and N=29 respectively 
are shown with full squares in Fig. ^. 



overestimate when ay is computed as a Gaussian quan- 
tity. 

Finally, the computed value for A accounts for an aver- 
age behavior of the velocity dispersion relative errors on a 
wide range of models and it can be considered as an upper 
limit for the relative errors, once given the number sample 
in bins. 

We have then compared the precisions on a single mea- 
surement, obtained as in Eq. (pl|), with those expected 
from the related distribution on 100 simulations and found 
that they are consistent within the errors. 

7.2. Errors on the rotation velocity 

We focus our analysis on the last spatial bin of the rota- 
tion velocity curves, because the outer regions are those 
for which we wish to derive the mass and angular mo- 
mentum estimate. The statistical errors for the velocity 
rotation measurements depend on the adopted model and 
bin dimensions. The bin dimension determines the num- 
ber statistics and involves the velocity range over which 
we compute the average rotation for a selected sample. 
The mass and rotation model determine the dispersion 
profile which is needed for the equilibrium, and this veloc- 
ity dispersion enters in to the error budget of the rotation 
measurement. In fact, using the definition of errors on the 
rotation velocity, we have 



SVrc 



1 



where ap — ayy^l + o-,nea/o'v — Cyay and 

SViot Cy ay 



Cy 

Vn 



(-) 

V (T / obs 



(22) 



(23) 
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Fig. 8. Precisions on rotation velocity. Left: behaviour of the generalised relative errors as function of and («/(T)obs, with 
the ±SD surfaces which are computed from the errors on the a and 6-parameters. Right: 15%, 20%, 30%, 40% contours of the 
generalised relative errors. 



Table 5. Here we list the set of parameter values for the 
fits displayed in Fig. |^ (right panel). 





Parameters in 


Sv/v fit 


iV 


A 


B 


29 


0.185 ±0.019 


0.00±0.03 


14 


0.27 ±0.05 


0.00±0.07 


13 


0.28 ±0.02 


0.02±0.03 


12.5 


0.29 ±0.02 


0.015±0.020 


11 


0.30 ±0.02 


0.05±0.03 


8.5 


0.38 ±0.03 


0.03±0.04 


8 


0.39 ±0.04 


0.02±0.05 


7 


0.50 ±0.04 


-0.02±0.06 



Here we want to check this dependence. For a fixed N 
(i.e. the mean number of particles in the last bin), one ex- 
pects the relative errors to be reproduced by the following 
function 



Af I \ -1 
— = A{v/a)ohs 

V 



B. 



(24) 



In Fig. ^ the best-fit for different N, as in Table | are 
shown. For these curves, B is always consistent with zero, 
while, from Eqs. (|2|) and (24), we expect A—Cy/Vn and 
Cv=1.08, for a typical crinca/o'v=0.4 in the last bins. To 
test this dependence, the A parameter obtained for the 
different error profiles in Fig. ^ is interpolated with the 
function a b: the best fit is for a=1.15±0.10 and 

b=-0.02±0.03, in agreement with the expected behavior. 

All these results were obtained for a wide range of mod- 
els, so they can be generalised and used to estimate the 
sample size needed to reach a targeted precision of the 
observed kinematical quantities. 

Once the kinematics is known from integrated light in 
the inner region of a galaxy system, we can estimate the 
minimum sample needed to reach a targeted precision in 
the outer regions, based on the inner (w/tTobs)- In Fig. ^ we 
show the dependence of rotation velocity precisions based 



on the results in Fig. |^, as function of (w/cr)obs and N. 
In the same Figure, the contours for different precisions 
are shown. The results of this analysis indicate that 1) 
the relative errors on the rotation velocity computed from 
a single measurement are consistent with those expected 
from the distribution of 100 simulation, 2) the precision 
on the rotational velocity is correlated with the internal 
kinematics of real systems and 3) this correlation is quite 
general and can be used in feasibility studies, once we need 
to estimate exposure times and plan observing proposals 
to reach a given sample size, i.e. a targeted precision. For 
example, if one needs to reach an accuracy on Vrot of 20%, 
at least 30 PNe are needed in the last bin for systems were 
{v/a > 1); to reach similar precision, larger samples are 
needed for systems with a lower value of {v/a). 

8. Summary and conclusions 

We have developed an algorithm to simulate real measure- 
ments of discrete radial velocity fields for a given mass 
model at equilibrium, using a simplified DF in a factor- 
ized form, and with different sample sizes. We use this 
algorithm to address two important questions when us- 
ing discrete radial velocity fields as a tool to study the 
kinematics of the outer regions of elliptical galaxies: 

1. Are the simple 3-parameter functional forms used to 
fit small sample of PNe radial velocities introducing 
biases on estimates of the kinematical quantities? 

2. How does the precision of a kinematical estimate de- 
pend on the sample sizes, taking into account the ro- 
tation and inclination of the system under study? 

The answers are the following: 

1) A bilinear fit to a sample coming from a system with a 
rotation curve like that of Centaurus A (Eq. ||) does intro- 
duce a bias in the velocity dispersion profile at large radii. 
This problem is important and pervades more sophisti- 
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cated non-parametric analyses. Any kind of smoothing al- 
gorithm applied to a small discrete sample of observed ra- 
dial velocities, as in the study of NGC 1316 by Arnaboldi 
et al. ( 1998D , is effectively fitting a more or less bilinear 
form to the data. Our simulations show that this is bound 
to introduce a bias in the estimate of the velocity disper- 
sion profile derived from the residual field, in particular 
for highly rotating systems. A flat rotation curve, or a 
rotation curve profile derived from averaged data in bins 
along the line of maximum gradient, should be adopted 
to derive the best estimate for a. Such a procedu re wa s 
indeed used in the case of NGC 5128 (Hui et al. |l995|) . 
The overestimate in a caused by fitting a plane can be up 
to 20% in a, i.e. 40% in the total mass. This bias corre- 
lates with the (u/cr) ratio, in the sense that higher rotation 
leads to stronger biases in a. We also identified the most 
reliably estimated quantity derived when using these sim- 
ple 3-parameter fields: the line of maximum gradient. This 
leads us to adopt the NFP, i.e. the analysis of the binned 
quantities along this P. A., as the best approach free of 
analysis-induced biases. 

2) We have found that the precisions obtained on single 
measurements are consistent with those obtained from the 
simulated distributions. Moreover we have found a gener- 
alised relation between the precision one can obtain for 
the observables in the outer regions as function of sample 
size and inner kinematics. This empirical relation can be 
very extensively used to plan and estimate observing times 
to study the dynamics of the outer haloes of giant early 
type galaxies with multi slit spectrographs like F0RS2 or 
VIMOS on VLT, or in slitless spectroscopy. 
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counts from the continuum galaxy background^ These 
quantities are related to the incoming fluxes as follows 



5007 



^csp Stci £ 



E, 



5007 



(26) 



is the flux received in the 5007 A line, t 



csp 



where i^5oo7 

is the exposure time, 5toi is the telescope surface, e is 
the total efficiency (instrumental efficiency -I- atmosphere 
absorption), i?5oo7 is the energy of the photons at A =5007 
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N, 
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A A tcsp •S'tci e 53(1(27) 



where Fy is the flux per unit area received by the con- 
tinuum galaxy background in the V-band, AA is the pass 
band of a typical narrow filter, and Ssd is the area of the 
seeing disk {Ssd = 7r(3 * FWHM/2.ib&f). 
Substituting ( p6| ) and ( ^7|) in the relation (^5|), we find 
that 



SNR = K 



5007 



^csp 5'tcl £ 



1 



£^5007 V -F5007 + -Fv AA5sd 
where 

^ -[^Q-0.4(8.3268(-;S-)'/'' + 12.7532-|-Aic-(B-V)) 

and 



5007 



_ '^Q-0A(vaii„-,+l3.7'i:) 



(28) 



(29) 



(30) 



(CiarduUo et al. 1989). In Fig. |, we plot the SNR 



for a 3.6m telescope, using Tcsp = 3 hrs, AA =60 A, 
FWHM=1.3" and e=0.5 and considering TOiim=27.2 mag, 
/^c=21.24 mag arcsec"^ and B-V=0.1 which are the typical 
values for an E galaxy in Virgo Cluster (Caon et al. 1994; 
McMillan et al. |1993D . The condition SNR{Riini) = 9 im- 
plies i?lim = 0.7i?o. 



Appendix: the theoretical signal-noise ratio (SNR) 

In our simulated velocity fields we introduce a cut-off in 
the projected PN density distribution determined by a 
inner limiting radius, at which the incompleteness of the 
observed PNe sample due to the bright continuum light 
in the central part of an E galaxy becomes significant. For 
a given luminosity profile, we consider all the PNe having 
magnitude m = muni which is the faintest magnitude in 
a PNLF for a complete sample, and we compute the SNR 
with respect to the galaxy background. We take a PN to be 
detect ed if it has a detected flux with SNR=9 (CiarduUo 
et al. 1987). Considering a negligible read-out-noise and 



sky background distribution, the SNR is defined as follows 



20 



0.5 



1 
R/a 



1.5 



Fig. 9. The R-coordinate dependence of the SNR. for our pa- 
rameter choice. 



SNR 



N 



5007 



^/NsoOT + A^gal 



(25) 



where A^5oo7 are the photon counts detected for a single 
PN in the 5007 A [0111] line and TVgai are the photon 



In Eq. (|2^), we are considering the inner bright part of 
galaxies where the photon counts from the star continuum is 
much larger than the sky background, i.e. A'gai >> Nsky In the 
same equation, Nsoor » Nsky by definition of detectability for 
a PN: SNR=9 implies that A^'soot > A''gai >> N^ky 
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Fig. 10. Velocity rotation curves for the systems with Vmax = O.SOiig. Top: Zl axis {edge-on); center: Z2 axis {edge-on); 
bottom: Zl axis (LOS: -45,-45,-45). The vertical long dashed line indicates iiiim, full points are the mean values obtained on 
100 simulations, solid line is the expected velocity rotation. Left: 500 PNe; center: 150 PNe; right: 50 PNe. 
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Fig. 11. Velocity rotation curves along Zl and Z2 for the systems with Fmax = O.SOvg, face-on case. Empty circles are the mean 
values obtained on 100 simulations along Zl, empty squares are the mean values along Z2, solid line is the expected velocity 
rotation. Left: 500 PNe; center: 150 PNe; right: 50 PNe. 
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Fig. 12. Velocity dispersion profiles for non-rotating systems: model without DM. Pull points are the mean values obtained 
on 100 simulations, solid line is the expected velocity dispersion profile. Top: 500 PNe;middle: 150 PNe; bottom: 50 PNe. 
Left: in radial bins; center: along Zl axis; right: along Z2 axis. 
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Fig. 13. Velocity dispersion profiles in radial bins: rotating model with DM and Vmax ~ 0.23wg, edge-on case. Full points are 
the mean values obtained on 100 simulations, solid line is the expected velocity dispersion profile along the Zl-axis, dot-dashed 
line is the expected velocity dispersion profile along the Z2-axis. The estimates by the fitting procedures (BP and PC) show the 
aaimuthal dependence while the NPP causes an overestimation of the velocity dispersion values. Top: 500 PNe; middle: 150 
PNe, bottom: 50 PNe. Left: bilinear fit procedure; center: no-fit procedure; right: flat-curve procedure. 
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Fig. 14. Velocity dispersion profiles for rotating model with DM along Zl (upper panel, Vmax ~ 0.30ug) and Z2 (lower panel, 
Knax = 0.23wg): edge-on case. Full points are the mean values obtained on 100 simulations, solid line is the expected velocity 
dispersion profile, empty squares are the estimates using dZmin- The bias introduced by the BF is indicated by the arrow (upper 
panel 500 PNe case; see also Fig. ^). In each panel: Top: 500 PNe; middle: 150 PNe;bottom: 50 PNe. Left: bilinear fit 
procedure; center: no-fit procedure; right: flat-curve procedure. 
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Fig. 15. Velocity dispersion profiles for rotating model with DM along Zl (upper panel, VInax ~ 0.30?;g) and Z2 (lower panel, 
Knax = 0.23115). LOS: -45:-45:-45. Symbols have the same meaning as in Fig. h4. The amount of the biases is negligible here 
because the (i'/cr)obs = Vrot/o"v ratio is lower. In each panel: Top: 500 PNe; middle: 150 PNe;bottom: 50 PNe. Left: bilinear 
fit procedure; center: no-fit procedure; riglit: flat-curve procedure. 
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Fig. 16. Velocity dispersion profiles in radial bins: rotating model with DM and Vmax = 0.23ug. LOS: -45:-45:-45. Symbols 
have the same meaning as in Figure |l^. Top: 500 PNe; bottom: 50 PNe. Left: bilinear fit procedure; center: no-fit procedure; 
right: flat-curve procedure. 
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Fig. 17. Velocity dispersion profiles in radial bins: rotating model with DM and Knax ~ 0.23i;g, face-on case. Here the expected 
velocity dispersion profile (solid line) has no azimuthal dependence and a perfect agreement is found with the estimates (filled 
squares). Top: 500 PNe;center: 150 PNe; bottom: 50 PNe. Left: bilinear fit procedure; center: no-fit procedure; right: 
flat-curve procedure. 



